Avoiding bias in reconstructing the largest observable scales from partial-sky data 
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Obscuration due to Galactic emission complicates the extraction of information from cosmological 
surveys, and requires some combination of the (typically imperfect) modeling and subtraction of 
foregrounds, or the removal of part of the sky. This particularly affects the extraction of information 
from the largest observable scales. Maximum- likelihood estimators for reconstructing the full-sky 
spherical harmonic coefficients from partial-sky maps have recently been shown to be susceptible 
to contamination from within the sky cut, arising due to the necessity to band-limit the data 
by smoothing prior to reconstruction. Using the WMAP 7-year data, we investigate modified 
implementations of such estimators which are robust to the leakage of contaminants from within 
masked regions. We provide a measure, based on the expected amplitude of residual foregrounds, for 
selecting the most appropriate estimator for the task at hand. We explain why the related quadratic 
maximum-likelihood estimator of the angular power spectrum does not suffer from smoothing- 
induced bias. 



INTRODUCTION 

It is unavoidable that we observe the Universe through 
the galaxy we inhabit. The foreground contamination in- 
jected by the Milky Way into full-sky cosmological data- 
sets must be modeled and removed, or the regions most 
conspicuously contaminated must be excised. Where no 
precise model of the foregrounds is available, cutting the 
sky is the most robust option, with the regrettable conse- 
quence that part of the signal is discarded along with the 
contamination. This includes information on the largest 
scales, which are valuable for a variety of reasons, includ- 
ing measurement of the integrated Sachs- Wolfe effect [l] 
and constraining primordial non-Gaussianity using trac- 
ers of large-scale structure [2] . 

It is impossible to uniquely recover the cosmological 
signal discarded in the sky cut. However, by writing 
down the likelihood for the region of the sky in which 
one trusts the data, it is possible to reconstruct an es- 
timate of the signal at large scales which maximizes the 
likelihood of the residual noise [3 . An alternative recon- 
struction scheme maximizes the posterior probability [4]- 
[8] of measuring the underlying cosmological signal given 
the available data and a prior theoretical expectation on 
the signal. 

The reconstructions estimate the large-scale (low-^) 
spherical harmonic coefficients, a^^, by treating the sig- 
nal at small scales as noise and only considering data 
external to the sky cut. If, as with the cosmic mi- 
crowave background (CMB), the field to be reconstructed 
is not band-limited, the proliferation of small-scale sig- 



nal makes the reconstruction noisy to the point of be- 
ing useless. Input maps are therefore smoothed - nec- 
essarily prior to cutting the sky - to truncate the signal 
and remove sources of confusion below a chosen angular 
scale [9 . However, smoothing leaks contamination from 
the masked region into the trusted data [lOj [11] , and the 
reconstructed spherical harmonic coefficients, airm are 
biased. In this work we explore the causes and expected 
magnitudes of this bias, and discuss how it can be miti- 
gated. 



MAXIMUM-LIKELIHOOD RECONSTRUCTION 

We begin with a description of the standard imple- 
mentation of maximum- likelihood CMB a^rn reconstruc- 
tion. The first step of the reconstruction process is to 
band-limit the temperature field by smoothing, typically 
with a Gaussian kernel of width 10° FWHM. As this 
removes information on the smallest scales, the map 
resolution can be downgraded to reduce computation 
time. The a^rnS in the range 2 < £ < fmax,rec (rep- 
resented for ease as the n^^^ -element vector a, where 
^atrr. = (^max,rec-l)(^max,rec+3)) are then reconstructed 
from the npix unmasked pixel temperatures, x, using [3] 



a = Wx. 



(1) 



The reconstructed spherical harmonic coefficients maxi- 
mize the likelihood of the residual noise, given the avail- 
able data, if the reconstruction matrix, W, is 
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W= [Y*C-^Y]-iY*C- 



(2) 



Here, Y are the spherical harmonics calculated at each 
unmasked pixeQ and C is the pixel-space noise covari- 
ance matrix 
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where R is uncorrelated, low-amplitude regularizing 
noise added to prevent C from becoming singular, Q is 
the smoothed theory CMB angular power spectrum, and 
Pi are the Legendre polynomials at unmasked pixels z, j. 



The sum over the multipoles £u 



< i < ^max ensures 



that the small-angular-scale CMB multipoles we do not 
wish to reconstruct are treated as noise. As stated above 
the CMB power must be artificially truncated to restrict 
the number of ambiguous modes accessible to the recon- 
struction. The smoothing kernel is deconvolved from the 
dimS after reconstruction by dividing the a^^s by the 
kernel's spherical harmonic transform. 

If a foreground signal bi is now introduced, so that 
X = Ya + b + n and a is the CMB signal uncorrelated 
with b and the noise n, the mean and variance of the 
reconstruction error e^rn = ^irn — Ciim are 



and 



(e) = Wb 



(ee')-(e)(e') = WCW' 



(4) 



(5) 



respectively. 

Throughout this work, we reconstruct the spherical 
harmonic coefficients up to ^max,rec = 10- The noise 
covariance matrix includes CMB power in the range 



< I < Ir. 



32 unless explicitly stated; this 



value is chosen such that modes with I > ^max are sup- 
pressed to (9 (few %) by the smoothing. The WMAP 
5-year best-fit C^s [12] are chosen for the theory CMB 
angular power spectrunj^ Input maps are smoothed at 
HEALPix [13] resolution A^side = 512 before being down- 
graded to Aside = 16 to retain the information required 
for the reconstruction while minimizing the number of 
pixels included in the noise covariance matrix. Diagonal 
regularizing noise R is added at the level of 2 /iK^ to al- 
low the inversion of the noise covariance matrix despite 
the presence of some null modes (which are irrelevant to 
the reconstruction). 



SMOOTHING-INDUCED BIAS 

This section outlines how a bias arises from smoothing- 
induced contamination of the unmasked pixels [TOj [11] . 



^ Without loss of generality, the reconstruction matrix in this work 

is formed from the real spherical harmonics. 
-^ Our results are not sensitive to the small differences between 

different WMAP releases in the best-fit cosmology. 



For clarity, we illustrate the smoothing-induced bias in 
the maximum-likelihood reconstruction (and, later, our 
proposed solutions) with results plotted in both harmonic 
and pixel space. 

It is beyond the scope of this work to estimate accurate 
foreground residuals resulting from different component 
separation methods; instead, we choose some residuals 
for the purpose of illustration. Following Ref. [14], the 
residual foregrounds are taken to be 1% of the differ- 
ence between the WMAP 7-year Internal Linear Com- 
bination (ILC) [15] and V-band [16] temperature maps. 
The resulting map is indicative of the extent and ampli- 
tude of the residual foregrounds in the WMAP 7-year 
ILC, the data-set for which the smoothing-induced bias 
was first described JOHTT]. It is important to note that a 
considerably higher level of contamination is present in 
the foreground-reduced maps for each individual WMAP 
frequency band provided by the WMAP team: contam- 
inants of ^ 50 times those used here are visible in these 
maps. The residual foregrounds are restricted to the pix- 
els within the sky cut, which in this work masks only the 
Galaxy and not individual point-sources. These degree- 
scale point-source cuts are, unsurprisingly, not found to 
significantly bias the large-scale reconstructed a^^s, and 
so for clarity Galaxy-only contaminants and masks are 
considered. 

The addition of simulated residuals and the spherical 
harmonic transform are linear, so the smoothing-induced 
bias is given by the reconstructed a^^s of the simu- 
lated residual map. We take the "standard" 10°-FWHM- 
Gaussian-smoothed maximum-likelihood reconstruction, 
using data outside the Galaxy-only part of the WMAP 7- 
year KQ85 mask, as our fiducial maximum-likelihood es- 
timator (hereafter the "Gaussian ML" a^rns). The Gaus- 
sian ML dirn^ generated from the simulated residual fore- 
grounds are plotted (deep-blue solid line) in Fig. [l] along 
with the full-sky airn^ (dotted line) for comparison. The 
smoothing-induced leakage from within the sky cut is 
clear to see: the reconstruction picks up about half of 
the power in the foreground residuals, even though the 
simulated residuals are entirely confined to the sky cut. 

The largest peaks in the bias affect the a^^nS satisfy- 
ing £ = 2n, 771 = for integer n (in Fig. IT] HEALPix 
index i'^ -\-i-\-m = {6, 20, 42 . . .}) [3 , and are positive for 
odd n, negative for even n. The pattern of these peaks 
can be explained by examining the reconstruction of the 
simulated Galactic residuals, plotted in Fig. [2| which 
are coldest along the Galactic plane. Smoothing these 
residuals reduces the pixel temperature values approx- 
imately symmetrically around the Galactic mask, and 
therefore pollutes the azimuthal modes which are also 
symmetric about the equator. The bias is positive for 
modes which have minima at the equator, and negative 
for those with maxima. There are also secondary peaks at 
^ = 2n + l,m = 1 (inFig.[l]^2^^ + m = {13,31,57...}), 
which again are positive for odd n, negative for even n. 
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FIG. 1. The spherical harmonic coefficients of the simulated 
foreground residuals, calculated using the Gaussian ML re- 
construction (deep-blue solid line), the full-sky data (dotted 
line) and the 10° Top-Hat ML reconstruction (light-blue solid 
line on x-axis). It is clear that the Gaussian ML reconstruc- 
tion leaks around half of the information from within the sky 
cut; this can be counteracted by smoothing with a top-hat 
kernel and using an extended mask. The HEALPix index 
i'^ -\-i-\-m maps each i, m combination to a unique index into 
the array of a^rnS. 
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FIG. 2. Clockwise, from top left: simulated Galactic fore- 
grounds, smoothed by a 10°-FWHM Gaussian and masked 
with the Galaxy-only KQ85 mask; the i < 10 Gaussian ML 
dimS reconstructed from the simulated Galactic foregrounds; 
the real spherical harmonics I4 and Ys 1 • The simulated fore- 
grounds yield a negative bias in the 140 mode, and a positive 
bias in the Ys 1 mode. 



These modes pick out the concentration of reconstructed 
foreground power in the Galactic centre. 

The leakage of information from within the sky cut can 
also be demonstrated in pixel- space. Taking the WMAP 
7-year ILC map, the full-sky a^^s are extracted, and the 
Gaussian ML reconstruction is performed. The spherical 
harmonic coefficients recovered in each case are then used 
to reconstruct the input ILC map using only 2 < ^ < 10, 
as plotted in Fig. [3J The maps formed from the full- 
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FIG. 3. The 7-year ILC 2 < i < 10 modes, plotted using 
(clockwise from top-left) full-sky a^mS, Gaussian ML a^mS, 5° 
Top-Hat ML a^mS, and Top-Hat WF a^mS. 



sky-a^^s (top-left) and reconstructed from Gaussian-ML- 
dim^ (top-right) are almost identical, even in the Galactic 
plane, confirming that the reconstruction has access to 
information well inside the sky cut. 

Although the bias illustrated in Fig. [l] looks problem- 
atic, its amplitude is at least reduced over using the 
contaminated full sky, and so the simple Gaussian-ML 
procedure may yet turn out to be useful. To compare 
with other possible approaches, we first need to dis- 
cuss its standard deviation - i.e. the scatter induced 
by the £ > ^max,rec modcs - which is calculated using 
Eq. |5] and plotted (as the deep-blue narrowest band) in 
Fig. |4] This is a few /iK at most. For a given -£, the 
modes that are reconstructed with the least precision are 
those with \m\ = ^, with the m = i modes typically 
the worst. This confirms the observation in Pontzen andl 
Peiris| ([2010^ [I4 that the sky cut removes the most infor- 
mation from modes with power concentrated towards the 
equator, and particularly those with extrema at = 0°, 
where the mask is at its widest [11 . The mask is plot- 
ted for reference in Fig. [5J along with examples of the 
affected modes. We see that, typically, three modes per 
£ will have increased bias or variance, but for most modes 
both the mean reconstruction error and its variance will 
be small. 

In order to compare estimators, we must first quantify 
their performance over the range of multipoles consid- 
ered. The performance measure 



= mean(e£^)^ + var(e£^ 



(6) 



provides the expected size of the reconstruction error e^^ 
for each mode: summing over all modes 
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(7) 



therefore yields a complete measure of each estimator's 
performance. Any alternative estimator which removes 
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FIG. 4. The one-standard-deviation ranges of the reconstruc- 
tion error e^rn for the Gaussian ML a^rnS (narrowest, deep- 
blue band), 5° Top-Hat ML a^mS (widest, hght-blue band), 
and Top-Hat WF a^mS (intermediate, mid-blue band). 
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FIG. 5. Top: the WMAP 7-year Galaxy-only KQ85 mask 
(light-blue central region) extended by 2.5° (mid-blue) and 
5° (dark-blue). Bottom: the real spherical harmonics 144 
and 14-4- The concentration of Yi±i mode power towards 
the equator results in increased estimator variance in those 
modes. 



the smoothing-induced bias should be preferred only 
if its Z value is lower than that of the Gaussian ML 
reconstruction, and, indeed, the contaminated full-sky 
airn^. In fact, using the 1% ILC-V foreground residu- 
als, Z c^ 33 /iK^ using the full-sky a^^s, compared to 
^ 265 /iK^ for the Gaussian ML estimator. If the residual 
foregrounds employed in this work are an accurate reflec- 
tion of those present in the WMAP 7- year ILC map, then 
the contaminated full- sky a^^nS provide a better estimate 
of the cosmological signal than the Gaussian ML recon- 
struction. The second form of Eq. [6] shows that both 
bias (see Fig. fl]) and variance (see Fig. H| in the recon- 
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FIG. 6. Top: the residual foreground map employed in this 
work (left), and a simple model comprising a 3°-wide — 20/iK 
equatorial band. Bottom: the full-sky a^rnS (dotted) and re- 
constructed Gaussian ML aimS (solid) of the residual fore- 
ground maps. The simple band model captures most of the 
features of the smoothing-induced bias injected by the more- 
complex residuals. 



struction increase the value of Z. The variance term is 
independent of the contamination, while the bias scales 
linearly with the contamination. Therefore, if the am- 
plitude of residual foregrounds in the ILC map is higher 
than in our illustrative example, the reduction in bias due 
to the use of the Gaussian ML a^rnS will eventually over- 
come the variance introduced by the reconstruction. For 
residual levels 3 — 4 times higher than those used here, 
the Gaussian ML reconstruction should be used instead 
of the contaminated full- sky a^^nS. However, as we have 
seen, the Gaussian ML estimator (as implemented thus 
far) does not eliminate the bias due to smoothing-related 
leakage of contaminants from within the masked region. 
While our simulated foreground residuals are simply 
meant to be indicative, we nevertheless expect that the 
smoothing bias is mainly sensitive to the amplitude of the 
residuals, and not their precise morphology. This can be 
seen in Fig. [6) where we have modeled the residuals as 
a simple bar in the Galactic plane, while rescaling the 
amplitude to match our 1% (ILC-V) model. This highly 
simplified model is able to capture most of the features 
of the bias in harmonic space, as seen in the lower panels 
of this Figure. 



ELIMINATING THE BIAS 

At this stage, we are presented with something of a 
conundrum: smoothing is essential to the reconstruction 
process, but it is exactly this smoothing that is biasing 
the results. The simplest solution to this issue is to re- 
move the areas of the sky that are within one smoothing 
scale of the main Galactic sky cut. However, the smooth- 



ing kernel typically used in the standard reconstruction 
algorithm is a Gaussian, with support across the full sky 
in pixel space, and the set of contaminated pixels is hence 
poorly defined. This problem can be solved by using a 
kernel with finite pixel-space support - for simplicity we 
choose a top-hat - as all contaminated pixels fall within 
a kernel's radius of the mask. 

Naively selecting the diameter of the top-hat smooth- 
ing kernel to be 10°, we reconstruct the simulated fore- 
ground residuals using the Galaxy-only KQ85 mask ex- 
tended by 5° - hereafter the "10° Top-Hat ML" recon- 
struction. The results are plotted in light blue (pale line 
along X-axis) in Fig. fl] the smoothing-induced bias has 
been eliminated. However, the measure of reconstruction 
quality has deteriorated dramatically to Z c^ 11 252 /iK^, 
significantly worse than the Gaussian ML reconstruction. 
There are two reasons for this increase in Z, which is now 
sourced entirely by increased variance in the reconstruc- 
tion. Firstly, the top-hat smoothing kernel has support 
over a greater range of multipoles than the Gaussian ker- 
nel, and so more ambiguous modes contribute to the co- 
variance matrix (in this case, and for all further kernels, 
we increase ^max to 1024 to capture all relevant modes, 
even though at 10° the smoothing kernel is effectively 
band-limited at £ '^ 100) |j The second factor is that the 
reconstruction-error variance increases very rapidly with 
the area of the sky that is masked [3] . This suggests the 
use of narrower smoothing kernels, although this neces- 
sarily increases the power of the high-^ noise. 

The interplay between the variance injected by de- 
creased smoothing and increased masking is shown in 
Fig. [7| Here, the total reconstruction-error variance 
is plotted for top-hat smoothing kernels of diameter 
3° — 10°, and hence mask extensions of 1.5° — 5°. The 
reduction in variance due to minimizing the sky cut dom- 
inates the added noise from narrower smoothing. We 
therefore select the width of the smoothing kernel to be 
as small as possible, given the resolution of the degraded 
map. At A^side = 16, the pixels are ^ 3° across, so 
to avoid injecting bias through pixelization [11] (which 
would not be captured by Eq. [sl we choose our optimal 
kernel diameter to be just larger: 5°. Hereafter, we refer 
to this reconstruction - using the 5° top-hat smoothing 
kernel and Galaxy-only KQ85 mask extended by 2.5° - 
as the "5° Top-Hat ML" reconstruction. 

The reconstruction-error variance calculated for the 5° 
Top-Hat reconstruction is plotted as the light-blue out- 
ermost region of Fig. [4J Even using the minimum pos- 
sible mask extension the reconstruction-error variance is 
still an order of magnitude larger than that of the Gaus- 
sian ML estimator. In terms of the measure of recon- 
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FIG. 7. The impact of mask extension and smoothing- kernel 
diameter on the quality of the bias-free maximum-likelihood 
reconstruction. Reconstruction-error variances are calculated 
using top-hat smoothing kernels of varying diameter, and 
extending the KQ85 mask by one kernel radius each time. 
Smoothing at the lowest-possible scale will produce the most 
faithful reconstruction. Note that the curve is not smooth 
as both the extent and shape of the mask change as it is ex- 
tended. 



struction quality, the biased Gaussian ML reconstruc- 
tion {Z ~ 265 /iK^) should be strongly preferred to the 
5° Top-Hat ML case (Z - 8466 /iK^) for residual lev- 
els comparable to those used in this work. As the Top- 
Hat ML reconstruction is unbiased, this value of Z is 
fixed (for a given sky cut). Thus, only if the residuals 
are greater than - 25% ILC-V wih the 5° Top-Hat ML 
reconstruction outperform the Gaussian ML reconstruc- 
tion. Note that the quality of the reconstruction could 
be improved further if it was performed at higher reso- 
lution, as smaller smoothing kernels could be used. This 
will necessarily have to be traded off against the increased 
computational requirements. 



REDUCING THE VARIANCE 

The increase in variance encountered when using ex- 
tended sky cuts is far beyond that expected due to the 
reduction in pixel count: reducing /sky from 81% to 74% 
should, assuming uncorrelated pixels for simplicity, in- 
crease the variance by only ~ 10%. The dominant issue 
is that the maximum-likelihood reconstruction allows the 
temperature field in the masked region of the sky to have 
infinite variance. For small sky cuts (and small ^max, rec) 
this is fine: one cannot "hide" large-scale power within 
the cut, and so the variance on the large-scale a^^s is low. 
Extending the KQ85 mask not only increases its overall 
width, but also closes a number of small gaps that allow 
the estimator limited access to the poorly-constrained 
equatorial modes (see Fig. p|. The estimator is therefore 



free to fill the cut with significant low--^ power (compare 
Fig. [3] top-right and bottom- left), and the estimator vari- 
ance rises rapidly. 

The variance of the reconstruction error can be reduced 
by enforcing a prior on the power within the Galactic 
cut using Wiener-filtering (Method 5 in Ref. [7] ; see also 
Refs. [IHSII^)- The Wiener-filtered a^rnS then maximize 
the posterior probability of reconstructing the underly- 
ing a^rnS, given the trusted data and theoretical power 
spectrum. In practice, this is achieved by adding a the- 
oretical (inverse) covariance matrix for the multipoles of 
interest to the reconstruction matrix 



W = [S-^ + Y^C-^Y]-^Y'C- 



(8) 



S here is a diagonal n^ 



X n„ 



matrix with elements 



equal to the smoothed theory power spectrum C^. The 
theory prior restricts this "Top-Hat WF" reconstruction 
from filling the sky cut with arbitrary power (see Fig. Is] 
bottom-right), and reduces the variance of the estima- 
tor's reconstruction error accordingly (see the mid-blue 
intermediate region in Fig. p|. The analytic expectation 
for the reconstruction-error variance is 



(ee*) -(€)(€*) = S[S + N]-iN, 



(9) 



where N is the noise power spectrum, defined to be 
N = [Y^C-^Y]-^ {i.e. the variance of the Top-Hat 
ML reconstruction error). By adding a theory prior to 
the power within the sky cut, and hence requiring finite 
power in that region, the Wiener-filtered reconstruction 
tends to produce a^rnS that are closer to zero than the 
maximum-likelihood case. While this could be seen as 
biasing the a^rnS toward lower valueq^ it can also be in- 
terpreted as being conservative, and applying the prior 
belief that the information within the mask is similar 
to the trusted information outside the mask. In other 
words, we should be happy to trade off a small multi- 
plicative bias against a significant reduction in variance. 
This is automatically encapsulated in the measure of 
reconstruction quality Z for the 5° Top-Hat WF estima- 
tor, which has improved to ~ 1521 /iK^. However this is 
still worse than that of the Gaussian ML estimator. The 
bias arising in the Top-Hat WF reconstruction is not from 
smoothing but from a prior, so Z is fixed for a given mask, 
and always lower than that of the corresponding Top-Hat 
ML reconstruction. For contamination levels of > 10% 
ILC— V (such as those found in the foreground-reduced 
maps for the individual WMAP frequency bands), the 
5° Top-Hat WF reconstruction therefore represents the 
most reliable estimator considered in this work. 



"^ Assuming for clarity zero noise, the Wiener-filtered reconstruc- 
tion yields a = S[S + N]~-'^a, ie. a multiplicative bias. Note 
that the ensemble average (a) = (a) = 0. 
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FIG. 8. The reconstructed WMAP 7-year ILC a^mS, calcu- 
lated using the Gaussian ML reconstruction (top) and the 
5° Top-Hat WF reconstruction (bottom). The shaded areas 
represent the estimator standard deviations. The modes that 
are most contaminated by the simulated foregrounds in the 
Gaussian ML reconstruction are indicated, along with their 
expected sign, by dashed (l^eveno) and dash-dotted (V^ddi) 
lines. 



We do not know the precise magnitude or morphology 
of the residual foregrounds in the WMAP data. We can, 
however, examine the Gaussian ML and Top-Hat WF 
reconstructions of the ILC by eye to determine if there 
is an obvious impact due to residual foregrounds. The 
dim^ generated from these two reconstructions are plot- 
ted iii^fjgj8J_alongwith the estimator standard deviation 
(\/(|^£mp) — {dim)'^)' Comparing the two plots, we see 
that there is very little difference between the dimS re- 
turned in each case. Further, no modes look statistically 
anomalous at the 3<j-level, even those that we expect to 
be contaminated from the simple residuals model used 
here. 



RELATION OF a^^ RECONSTRUCTION TO THE 
QML ESTIMATOR FOR THE C^S 

We have so far discussed estimating the fuh-sky a^rnS 
from cut-sky data, which is equivalent to reconstructing 
the smoothed temperature field. However, in the context 
under which the smoothing-induced bias was revealed [9]- 
[llj it is in fact only the angular power spectrum Ci of 
the temperature field which is required. 

A popular method for estimating the full-sky angu- 
lar power spectrum is to adopt the quadratic maximum- 
likelihood estimator as first derived in Ref. [17]. It has 
been noted (see e.g. Section 3 of Ref. [9] for a complete 
discussion) that the QML estimator can be formed using 
the maximum-likelihood airn^. On the surface, the QML 
estimates (henceforth denoted C^) may therefore seem 
to be susceptible to similarly problematic contamination 
from a smoothing stage. 

However, this is not the case: in fact the C^s are far 
more robust to the content of the cut because the smooth- 
ing can be conducted on vastly smaller scales (e.g. 1° in 
Pontzen and Peiris| ( |2Q1QD [Ti]). Note that |Copi et c 
(2Qll)[ll| miss this point, because they consider only 
two extreme cases: (i) smoothing at 10° and (ii) failing 
to smooth. They therefore reach the erroneous conclu- 
sion that the QML estimator is susceptible to contami- 
nation from within the mask. We explicitly verified that 
the pipeline used by [Pontzen and Peiris| ( [2QlQ| ) [14] is in- 
dependent of any contamination placed fully inside the 
mask. 

The above paragraphs at first appear to be contra- 
dictory, since they simultaneously claim (a) that the 
QML power spectrum estimates can be formed out of 
the ML temperature field reconstruction; and (b) that 
the QML power spectrum estimates can still be con- 
structed from maps smoothed on degree scales (whereas 
the a^rnS will necessarily become noisy for sufficiently 
high ^max,rec)- Howcvcr, this is not a true contradic- 
tion because the QML estimates are not formed directly 
from the noisy a^^s, but rather through an expression 
(Eq. 23 of Ref. ^J which specifically downweights poorly 
constrained modes. It is this cautious treatment of am- 
biguous modes which makes power spectrum estimation, 
as opposed to a£rn reconstruction, so well-behaved, irre- 
spective of the shape of the smoothing kernel employed. 



DISCUSSION 

Maximum-likelihood estimators, dirn^ are often used 
to reconstruct the large-scale spherical harmonic coeffi- 
cients, a£rri^ from partial-sky data. The technique relies 
on smoothing to restrict the amount of small-scale noise 
accessible to the reconstruction, but smoothing has been 
shown to contaminate "clean" pixels with residual fore- 
grounds from within the sky cut. In this work, we have 



examined the impact of this smoothing-induced bias on 
the maximum-likelihood reconstruction. We have shown 
that it is possible to mitigate the bias by removing the 
contaminated regions, but these are only well-defined if 
smoothing is performed using a kernel with finite sup- 
port on the sky. This precludes the use of the commonly 
used Gaussian kernel. Cutting a larger portion of the sky 
greatly increases the variance of the reconstruction, but 
it is possible to counteract this effect by enforcing a prior 
on the reconstructed coefficients using a Wiener filter. 
We have therefore proposed an estimator - using top-hat 
smoothing, extended masks and a Wiener-filtered recon- 
struction - which does not suffer from smoothing-induced 
bias. By considering the expectation of the square of the 
reconstruction error, Z = ^^ ^{{ciim — ciimY)-, we have 
compared the performance of the maximum-likelihood 
and Wiener-filtered estimators in the presence of sim- 
ulated CMB foreground residuals. 

The reconstruction performance measure Z scales 
with the estimators' bias and variance, which in turn 
are governed by the amplitude of contamination and 
the size of the sky cut, respectively. The fiducial 
maximum-likelihood reconstruction is performed using 
relatively small sky cuts, but is susceptible to con- 
tamination through smoothing-induced bias; the finite- 
smoothing Wiener-filtered reconstruction does not suffer 
from smoothing-induced bias, but makes use of extended 
masks. Increasing the level of contamination therefore 
increases Z for the maximum-likelihood reconstruction 
only, which suggests that there is a level of contamina- 
tion above which one should switch from the maximum- 
likelihood to the Wiener-filtered reconstruction. 

Given an estimate of the morphology and amplitude 
of the contaminants within the cut sky, one can predict 
which modes will be biased and by how much, and hence 
determine the threshold at which one should swap esti- 
mators. We find that this threshold is relatively insensi- 
tive to the precise morphology of foreground residuals at 
large scales, and is mainly governed by their amplitude. 
Calculating Z for the two estimators in the presence of 
estimated foreground residuals, we determine this thresh- 
old to be ~ 10 times the amplitude of the foreground 
residuals used in this work. Assuming that the ILC con- 
tains similar levels of contamination to those used here, 
we therefore recommend the use of either the contami- 
nated full- sky airaS or the fiducial maximum-likelihood 
airn^ when handling this data-set. However, when us- 
ing foreground-reduced maps for individual WMAP fre- 
quencies, which contain much greater foreground resid- 
uals, the Wiener-filtered reconstruction will provide the 
best estimate of the large-scale underlying CMB signal. 
Note that, as the Wiener-filtered a^^nS are a maximum- 
posterior solution, care must be taken if the reconstruc- 
tion output is being used for further model-selection 
steps. The reconstruction techniques are, however, most 
commonly used to test the null hypothesis, in which case 



the prior employed in this work is completely appropri- 
ate. 

For problems requiring only a power spectrum (as op- 
posed to the full temperature field) the issues described in 
this paper are essentially irrelevant because the smooth- 
ing can be conducted on vastly smaller scales, the result- 
ing range of poorly constrained modes being automati- 
cally downweighted. 
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